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We use the WMAP 7 data to place constraints on oscillations supplementing an almost scale- 
invariant primordial power spectrum. Such oscillations are predicted by a variety of models, some 
of which amount to assuming there is some non-trivial choice of the vacuum state at the onset of 
inflation. In this paper we will explore data-driven constraints on two distinct models of initial 
state modifications. In both models the frequency, phase and amplitude are degrees of freedom of 
the theory for which the theoretical bounds are rather weak: both the amplitude and frequency 
have allowed values ranging over several orders of magnitude. This requires many computationally 
expensive evaluations of the model CMB spectra and their goodness-of-fit, even in a Markov Chain 
Monte Carlo (MCMC), normally the most efficient fitting method for such a problem. To search 
more efficiently we first run a densely spaced grid, with only 3 varying parameters; the frequency, the 
amplitude and the baryon density. We obtain the optimal frequency and run an MCMC at the best 
fit frequency, randomly varying all other relevant parameters. To reduce the computational time 
of each power spectrum computation, we adjust both comoving momentum integration and spline 
interpolation (in I) as a function of frequency and amplitude of the primordial power spectrum. 
Applying this to the WMAP 7 data allows us to improve existing constraints on the presence of 
oscillations. We confirm earlier findings that certain frequencies can improve the fitting over a model 
without oscillations. For those frequencies we compute the posterior probability, allowing us to put 
some constraints on the primordial parameter space of both models. 



I. INTRODUCTION 

The observed statistical distribution of temperature 
fluctuations in the cosmic microwave background (CMB) 
is believed to be largely determined by the physics in 
the very early Universe. These CMB fluctuations were 
sourced by quantum fluctuations during an epoch of ac- 
celerated expansion early on in the history of the Uni- 
verse, known as inflation. They then induce curvature 
perturbations in the geometry of spacetime, which are 
preserved after horizon crossing during inflation. Once 
they re-enter the horizon at some later time, they cou- 
ple to radiation and matter, becoming responsible for 
the observed statistical distribution of the CMB tem- 
perature fluctuations and the large scale structure (LSS) 
of matter in the Universe, respectively. Within the 6- 
parameter ACDM these initial conditions are described 
by only 2 parameters: the amplitude of the primordial 
power spectrum of scalar perturbations, and the tilt n s 
describing the scale-dependence of the power spectrum 
to first order. By investigating the CMB power spec- 
trum, one is therefore able to probe high energy physics 
in the early Universe. Given the large number of mod- 
els describing the physics of the early Universe (see e.g., 
lj for a recent overview), these two parameters (within 
6-parameter ACDM) are not sufficient to be able to dis- 
criminate between various proposed models. Additional 
degrees of freedom, derived from the statistical analysis 
of the late time distribution of temperature (CMB) or 



density (LSS) fluctuations, could potentially be used to 
break the degeneracy between various models. Possible 
extensions to the 6-parameter model include tensor de- 
grees of freedom (gravitational waves), higher order cor- 
rections to the scalar and tensor power spectra, and de- 
viations from Gaussianity measured through high order 
correlation functions. The ultimate goal, of course, is to 
formulate a theoretical prediction of what these degrees 
of freedom should be, and constrain these from the ob- 
servational data, just as we constrain the 2 parameters 
describing the tilt and the amplitude of the primordial 
power spectrum of scalar perturbations. The additional 
constrained parameters test the uniqueness of a proposed 
model and thus contribute to the understanding of the 
physics governing the early, inflationary Universe. 

In this paper we will consider modifications to the pri- 
mordial power spectrum. In particular, we will search for 
evidence of oscillations in the almost scale-invariant spec- 
trum. The motivation to search for these modifications 
is provided by a rapidly increasing number of theoretical 
models expecting such features in the primordial corre- 
lation statistics: e.g., in the power spectrum [2"rH2l [Til- 
[151 |2"2"] and bispectrum [T?r[2"2"] . In this paper we will 
only consider oscillations in the power spectrum. |23j 
have attempted to constrain features in the bispectrum 
and [23] have proposed a method to effectively search for 
features in the CMB bispectrum; recently LSS data has 
been proposed to search for these type of bispectra [25] . 
The frequency, amplitude, and possibly the phase and 
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first order scale dependence of these features are deter- 
mined by the detailed physics of inflation. Detecting or 
constraining these parameters would help us determine 
the precise physics of inflation. 

This is not the first time features in the power spec- 
trum will be explored. Notably 0[26], [23 EH], [29], [30], 
[3"T] . and []~2] have all investigated possible oscillations 
in the primordial power spectrum. Although there has 
been no clear detection, the data certainly allow small 
oscillations, as can be seen from fig. [T] Previous analyses 
were done on WMAP 1-, 3-, and 5-year data (as well as 
SDSS data). In this paper we will aim at extending and 
improving the analysis using the latest WMAP release, 
WMAP 7 [22] • One of the key hurdles in searching for os- 
cillations in the data is the frequency of the oscillations. 
The frequency of the oscillations in the primordial power 
spectrum is a free parameter which spans several orders 
of magnitude for many of the proposed models. When 
probing the joint likelihood of our cosmological model, 
the large range in frequencies results in a number of is- 
sues [33] , 

Foremost, it requires a high resolution in 'sample 
space', i.e., looking for the best fit parameter to the 
data requires a large number of computations of the 
CMB power spectrum from the primordial power spec- 
trum (using, e.g., CAMB [3J] or CMBfast [35]). This ne- 
cessitates an efficient computing scheme. Unfortunately, 
the computation of the late time power spectrum from 
the primordial one involves a convolution of the trans- 
fer function A;(fc) with the primordial power spectrum 
P(k). The appearance of rapid oscillations in both the 
transfer functions as well as the primordial spectrum re- 
quires smaller and more frequent steps in the integration 
in comoving momentum space fc, increasing the required 
computational time for each run significantly. 

Secondly, the late time spectrum must to be computed 
for each angular scale I. Usually, given the smoothness 
of the primordial spectrum, it suffices to do a spline in- 
terpolation on a number of Ci. This reduces the com- 
putational time, since the computation of the transfer 
functions 1 is the most time-consuming part. With the 
addition of oscillations on top of the smooth primordial 
spectrum, the number of C\ necessary to obtain an ac- 
curate fit of the interpolated C/ will depend on the fre- 
quency of the primordial oscillations. As this frequency 
increases, at some point all I will need to be considered 
in order to resolve the superimposed oscillations. Com- 
puting all I requires us to compute all transfer functions, 
which again increases the computational time. Lastly, 
there will often be a number of frequencies that tend to 
improve the fit within the large range we explore, rather 
than just a single one. For example, if a frequency uj is 



a good fit to the data, there is a fair chance 2cu will be 
a good fit as well. This is a major issue, as a multidi- 
mensional parameter space is most effectively searched 
through a Markov Chain Monte Carlo (MCMC), which 
is a random process. Therefore, when the frequency is 
not fixed, an MCMC approach to the best fit is not effi- 
cient as the likelihood is spiked and the random nature of 
the MCMC chain will lead to frequent tunneling of one 
local maximum to another within the multidimensional 
joint likelihood. |12j showed a way to circumvent this 
issue by first taking a large number of samples fixed on 
a grid. Obviously a grid does not allow us to vary all 
parameters within our cosmological model, as the num- 
ber of samples grows quickly with the dimensionality of 
parameter space. Instead, the priority lies in varying the 
initial conditions, i.e., the oscillatory component of the 
primordial power spectrum. Once the best fit has been 
determined (within the prior frequency range, and the 
grid resolution), one can run MCMCs with a fixed best 
fit frequency determined through the grid. The joint like- 
lihood is expected to no longer contain local maxima and 
an MCMC should converge quickly. This eventually al- 
lows us to put constraints on the amplitude and perhaps 
the first-order scale dependence of the amplitude of the 
oscillatory feature. 

Although there are a large number of different features 
and oscillations predicted by a variety of models, in this 
paper we will focus on only two, distinct, theoretically 
motivated modified primordial power spectra. We will 
introduce these modified power spectra in section In] In 



section III we explain how to optimize the search by mak- 
ing the numerical computation of the power spectrum 
frequency-dependent. To find a best fit frequency before 
we apply the MCMC, we use grid sampling. We report 



1 These can generally not be precomputed, as they depend on all 
relevant parameters of the theory, which are varied when search- 
ing for the best fit. 



our findings in section IV Once we have established a 
best fit value of the frequency we run an MCMC with 
that best fit value for two models. The results are dis- 
cussed in section[V]and we compare these to theoretically 
derived constraints on primordial parameter space. We 
conclude in section IVTT1 



II. TWO MODELS OF INITIAL STATE 
MODIFICATIONS 



Although the standard BD vacuum state during infla- 
tion is an excellent fit to the currently available CMB 
data, theoretical considerations have questioned its va- 
lidity and uniqueness over the last decade (see [2HS] and 
references therein). The main reason for casting doubt 
on the BD assumption is the fact that the temperature 
fluctuations in the CMB ultimately find their origin in 
quantum fluctuations in a vacuum state right at the onset 
of inflation. Predicting the spectrum of CMB tempera- 
ture fluctuations requires determining this initial vacuum 
state for quantum modes at extraordinarily high mo- 
mentum scales, far beyond any high-energy cut-off scale, 
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where a perturbative quantum field theory description 
is no longer expected to be accurate. From a theoret- 
ical point of view it is not understood why the time- 
dependent inflationary background would transfer this 
unknown high-energy physics to a regular quantum field 
theory description in the BD state at later times. In wait 
for a more complete understanding this vacuum state 
ambiguity has motivated phenomenological approaches, 
usually relying on low-energy effective field theory expec- 
tations, that have suggested the presence of high-energy 
corrections to the BD state. A typical prediction of these 
models is the appearance of (small) oscillations on top of 
the standard primordial power spectrum of inflationary 
perturbations. 

Besides constraining these characteristic oscillations in 
general, an additional goal of this work is to study to 
what extent the currently available CMB power spec- 
trum data can distinguish between two rather general 
classes of models that have been proposed to describe 
initial state modifications. The first class of models is 
known as the Boundary Effective Field Theory (BEFT) 
approach to determine initial-state modifications [7, Ell- 
in this proposal one fixes an initial time where one calcu- 
lates corrections to the usual BD initial condition using 
a low-energy effective boundary Lagrangian. The result 
is an explicitly scale-dependent Bogolyubov parameter 
/3k describing the modification with respect to the usual 
BD vacuum state, giving rise to strongly scale dependent 
oscillatory corrections to the primordial power spectrum. 

The second class of models is known as the New 
Physics Hypersurface (NPH) approach to initial-state 
modifications [36] . In the NPH scenario one traces 
every momentum mode back to some large physical scale 
of new physics M and imposes, on a rather ad-hoc basis, 
the standard flat space vacuum state (corresponding to 
positive frequency modes only), mode by mode, resulting 
in a k independent Bogolyubov parameter /3k- A small 
departure from scale-invariance only arises after taking 
into account the slow-roll evolution of the Hubble pa- 
rameter, which affects the amplitude of the oscillatory 
corrections described by /3k oc H/M. A recent effort 
by Jackson and Schalm 16J to compute the low-energy 
vacuum state effects of a massive scalar field in an in- 
flationary background carries the important prospect to 
ground the NPH proposal in a more solid effective field 
theory description. The authors showed that integrating 
out an arbitrary massive scalar field can indeed affect 
the vacuum, resulting in rather similar corrections to the 
inflationary power spectrum as in the original NPH sce- 
nario. For our purposes the power spectrum predicted 
in [16] will be used to set-up a general phenomenological 
parameterization of the almost scale-invariant NPH class 
of models. 

These two classes of models describing initial state 
modifications, NPH and BEFT, serve as a nice bench- 
mark to study how constraining and distinguishing the 
most recent CMB power spectrum data is when consid- 



ering superimposed oscillations. The data analysis is in 
fact best performed using a primordial power spectrum 
parameterization of the oscillating corrections involving 
a set of independent parameters. The general parame- 
terization of the primordial power spectrum for the two 
classes of models that we will make use of are [5] 

P(k) = P (fc) [1 + f3kA Q sin (2A 1 k + (/))} (1) 

for the BEFT scenario, while for the new Jackson and 
Schalm scenario [TH] (from here on NPH) 

P(k) = Pq [o + a x lnfc/fc*+ 

[ol2 + 03 \nk/k*) sin (04 lnfe/fc* + £)] . (2) 

In eq. Q /3 is a parameter determined by the details of 
the physics in the ultraviolet, and is naturally expected to 
be 0(1), while A a and A\ are (related) parameters pre- 
dicted by the BEFT method. Po(k) is the power spec- 
trum from canonical slow-roll inflation in a BD state. 
This contains a possible tilt n s . In eq. pi, P * is the 
scale independent power spectrum, and a possible tilt is 
not included. The reason for this difference is that the 
oscillatory correction in the NPH scenario can have a tilt 
03 that differs from the tilt a\ = n s — l. For the BEFT we 
will be interested in possibly constraining Aq, A\ and </>, 
while for the NPH model we will consider 02, 03, 04 and £. 
The parameters a and ai will be considered zero-order 
contributions to the power spectrum and are constrained 
through n s and P*. As already briefly alluded to, we 
will constrain each parameter independently and apply 
theoretical dependencies between those parameters only 
afterwards. For example, A\ in the BEFT model is re- 
lated to the amplitude Aq. The reason not to implement 
this directly is because separating these parameters re- 
sults in a much 'cleaner' sampling, because the frequency 
causes the likelihood to be highly irregular. We will apply 
the theoretically predicted relations between the ampli- 
tude and the frequency once the posterior distribution 
has been determined. 

For more specific details and theoretical background on 
the NPH and BEFT models we refer to [9] and QI] , 
respectively. 

III. CODE OPTIMIZATION 

We have modified the CAMB code [34] in order to ef- 
ficiently search for oscillations. The modification is built 
upon work done in [571 HE] • The late time power spec- 
trum basically is a convolution of the transfer function 
A;(fc) and the primordial power spectrum P(k) 

/>oo 

Ci oc k 2 dkP(k)Af(k). (3) 
Jo 

The transfer function contains all the physics that gov- 
erns the evolution of the Universe and describes an in- 
terplay between radiation pressure and gravitational col- 
lapse. This causes the transfer function to be highly oscil- 
latory. Consequently, to gain sufficient accuracy one has 
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FIG. 1: WMAP 7 data with added errors from measurement 
and cosmic variance. The error bars are derived from the diag- 
onal terms in the Fisher matrix. The multipole moments are 
slightly coupled (inducing off-diagonal elements) , so a correct 
treatment of errors requires use of the entire Fisher matrix, 
which is done when calling the likelihood WMAP code. 

to numerically integrate over fc-space with an adequate 
number of steps. The transfer functions are not only a 
function of comoving momentum, but also depend on the 
angular scale I. To increase the speed of the code, one 
usually does not compute the complete transfer matrix, 
but rather limited number of I, and spline interpolates 
between them. 

If the primordial power spectrum is smooth both the 
integral over comoving momentum and the spline inter- 
polation in I can be done with limited samples. Once 
we allow for (rapid) oscillations in the primordial power 
spectrum, the number of samples needs to be increased. 

Numerically, the most time consuming application is 
the number of transfer functions we will have to com- 
pute. Therefore the fewer transfer functions we compute, 
the faster we can determine the angular power spectrum 
Cj. Obviously, if the number of oscillations in our pri- 
mordial spectrum is large, we will need many transfer 
functions to resolve those using a spline interpolation. 
If the wavelength of the primordial signal is 61 ~ 1, no 
spline interpolation will resolve the primordial signal. 



A. Frequency-Dependent Power Spectrum 
Computation 

We have modified CAMB and COSMOmc in order 
to optimally compute the power spectrum given a pri- 
mordial frequency. Fixing all parameters to the best fit 
WMAP 7 values [32], except for the frequency in the 
primordial power spectrum, we investigated the conver- 
gence in the late time power spectrum when we increased 
the sampling in fc- and Z-space. As a null hypothesis, 
with 'optimal' accuracy, we took the primordial spectrum 
without oscillations. By plotting the frequency against 



FIG. 2: Sampling increase for splining in multipole space as 
a function of Ai for the BEFT model. We fixed A = 10, 
its largest possible value. As the frequency A\ increases, it 
requires an increasing number of neighboring I in the inter- 
polation to maintain the same accuracy (which was set by 
WMAP best-fit standards). Around Ai ~ 2000 the sampling 
factor increases beyond 50, which implies that all I's are re- 
quired (i.e., splining is not sufficient and one needs to compute 
the full transfer matrix). 



the sample increase required to retain the same accu- 
racy, we obtained an estimate of the optimal sampling 
for every frequency. As an example, in fig. [2] we show 
the increase in /-sampling required to obtain an accurate 
computation of the power spectrum as a function of the 
BEFT frequency A±. Only for relatively low frequencies 
A\ < 2000 and 04 < 60 can we use splines in Z-space to 
determine the CMB power spectrum. At higher frequen- 
cies we need to compute the full transfer matrix in order 
to resolve oscillations 2 . For the sampling in fc-space the 
optimal sampling was derived when computing all trans- 
fer functions. We find that over the range of frequencies 
we will probe in this paper we need to approximately 
double the resolution in fc-space when doing the integra- 
tion. 



IV. GRID SAMPLING 

Given the number of free parameters describing the ini- 
tial conditions (at least A s , n s and the frequency and am- 
plitude of the oscillating correction) , many computations 
of the power spectrum are required to obtain the best fit 
(distribution). A commonly used approach is to apply 



2 The sampling is from a set of predetermined I (e.g., I = 
{1, 2, 3, ....100, 150, ...1000, 1100, ...}). The number of points used 
in the interpolation is increased within this set of predetermined 
I until this set is unable to achieve the same accuracy as one 
would obtain by computing all I. This set is based on efficient 
computation of a smooth primordial spectrum. 
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FIG. 3: (NPH) The 68% and 95% confidence levels for amplitude versus the frequency 04 for the low frequency grid 
1 < 04 < 200 marginalized over the baryon density. There are many local peaks in the marginalized likelihood. The two most 
likely grid points are (04 = 46.5, 02 = 0.056) and (04 = 98.41, a2 = 0.147). Towards higher frequencies, larger amplitudes are 
allowed by the data but do not necessarily represent significant improvements of the overall fit. 




FIG. 4: (NPH) The 68% and 95% confidence levels for a part of the high frequency grid. Note that the absence of a match 
between the low and high frequency regimes is artificial. In the high frequency grid, there are only 6 samples for 100 < 0,4 < 400. 
Additionally, a 100% likelihood was assumed within one grid. For very high frequencies, the wavelength of the oscillations is 
too small to be resolved at large angular scale (I < 200). Primordial oscillations are resolved only in the second peak and 
beyond. Very high frequencies can produce glitches in the large scale wing of the first peak in the angular power spectrum. For 
some specific frequencies and amplitudes these glitches can match outliers in the wing of the first peak of the observed angular 
power spectrum resulting in A\ 2 — 12 for 04 ~ 998. 



a Markov Chain Monte Carlo (MCMC), randomly pick 
values of a set of parameters, and based on the likelihood 
of the computed run, reject or accept this point to be part 
of our parameter probability estimate. This approach is 
highly efficient, as we only compute a limited number of 



points in the multidimensional likelihood to determine 
the posterior distribution. Once we add oscillations to 
the primordial power spectrum an MCMC method be- 
comes less efficient, as the likelihood is expected to be- 
come irregular in the coordinate of the oscillation (the 
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frequency) and there are many local maxima in the like- 
lihood function. At some point, we will leave this local 
maximum and end up in another. Therefore constraints 
on the initial conditions are hard to recover since the 
MCMC will constantly move to different local maxima 
in likelihood space 3 . proposed a different approach. 
Instead of an MCMC over all parameter space, they first 
considered a grid over a limited number of parameters. 
The parameters to vary in the grid should be those that 
determine the primordial conditions (e.g., the frequency). 
They also identified the baryon density Of, to be degener- 
ate with the amplitude and frequency of the oscillations, 
as it influences the height of the first peak. Finally they 
showed that oscillations in axion monodromy models of 
inflation are a better fit primarily to the first peak in the 
angular power spectrum, and hence could mimic some 
of the effects produced by the baryon energy density fif,. 
The grid therefore samples at least 3 parameters (am- 
plitude, frequency and baryon density). The other pa- 
rameters are fixed to the WMAP 7 |32] best-fit values. 
As such the null hypothesis (no oscillations) is the fit to 
beat. The advantage of the grid is that one probes the 
likelihood completely, although under fixed conditions for 
most of the ACDM parameters. Under the assumption 
that these other parameters are not (strongly) correlated 
with the varying parameters, we should be able to de- 
termine the absolute maximum in the likelihood rather 
than a local one. Once this local maximum has been de- 
termined, we can perform an MCMC with the frequency 
set to the best fit value, and allow all other parameters 
to vary. This should probe the likelihood of, e.g., the 
amplitude more efficiently. 

We ran three grids, one for the BEFT model and two 
for the NPH model. While in this paper we have opted 
to avoid details about the exact theoretical prediction of 
the parameters in each model, we decided to consider 
two grids for the NPH model focusing on two frequency 
regimes 04. This is motivated by the fact that the NPH 
model is not expected to have a significant number of os- 
cillations, as the theoretically predicted frequency is slow- 
roll suppressed, i.e., 04 cx e. We have therefore investi- 
gated low frequencies (f < a 4 < 200) with 200 log-spaced 
samples, and a high frequency regime (100 < 04 < 1000) 
separated into 500 logarithmically spaced samples. The 
amplitude 02 runs from to 0.4 in 120 and 200 equidis- 
tant steps, respectively. For the NPH grid = ( = 0. 
For the low frequency grid we set 0.021 < ^h 2 < 0.026 
in 10 equidistant steps, resulting in a total of 240.000 grid 
points, while for the high frequency regime we considered 
0.02 < Q,iyh 2 < 0.027 in 16 equidistant steps, with a total 
of 1.6 million grid points. 

In figures [3] and [4] we show the confidence contours for 
grids obtained for the NPH model of the frequency ver- 



Constraints on other parameters unrelated to the initial condi- 
tions can be recovered with sufficiently large samples. 



sus the amplitude in the low and high frequency regime, 
marginalized over the baryon energy density Peaks 
are areas in which the fit is best for non zero values of 
the amplitude of the modification, while valleys represent 
frequencies which are not a good fit to the data and the 
best fit is no modification. For example one can consider 
a likely frequency (peak) and plot the probability of the 
amplitude for that frequency to find that the most likely 
value for the amplitude is non zero. In fig. [5] we show the 
joint likelihood contour plot for the amplitude 02 and 
the baryon density Q.^h? marginalized over the frequency 
04. We find that 02 = is the most likely value in the 
low frequency grid. We also derived the joint likelihood 
for one of the best-fit frequencies (a 4 = 98.74) to show 
that the best-fit point has a non-zero amplitude ai with 
almost 95% CL. 

For the high frequency regime we show the effect of the 
highest frequencies on the marginalized amplitude a-i in 
fig. [6] Note that the confidence levels are determined 
assuming that the grid contains all possible values the 
frequency could have. In other words, one would hope 
that for either extremely large or small frequencies the 
likelihood of the fit would go to zero. The problem is 
that it does not, and therefore we can only determine 
the confidence levels within a prior determined parame- 
ter domain. We have partly motivated this domain on 
theoretical arguments. Observationally, data is the lim- 
iting factor. 

It should be obvious that many different frequencies 
represent good fits. As we had foreseen, this complicates 
running a large MCMC for all relevant cosmological pa- 
rameters. The best fit point we find for the NPH model 
(A% 2 12) centers around very high frequencies, with 
04 = 980 and a-2 — 0.39. In fact for this frequency, the 
best fit amplitude probably lies outside the domain of 
< a2 < 0.4. This implies a relatively large number of 
primordial oscillations. The angular power spectrum for 
such high frequencies has most of its oscillations damped, 
since this frequency is only resolved at scales beyond the 
first peak. [27] showed that for these scales the ampli- 
tude will be suppressed, similar to how the overall power 
is damped. Given the large measurement error at these 
scales, we would not expect the fit to improve that much. 
It turns out that at large angular scales, the barely re- 
solved high frequency causes glitches in the large scale 
wing of the first peak. The WMAP data contains out- 
liers in the large scale (small I) wing of the first peak 
(attempts have been made to understand this in terms 
of a feature in the primordial spectrum, see e.g. [29] and 
recently a principal component analysis has been per- 
formed to search for such a feature in [? ]); these are 
fitted for very specific unresolved frequencies of the pri- 
mordial power spectrum. These glitches are expected as 
the angular power spectrum will no longer be able to 
sample all oscillations appearing in the primordial power 
spectrum. We have not been able to get rid of these 
glitches by increasing the overall accuracy of the numer- 
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FIG. 5: (NPH) The 68% and 95% confidence levels of tlbh 2 versus a^. The left shows the joint likelihood after marginalizing 
over the frequency 014 (low frequency grid), while the right shows the joint likelihood for the best fit frequency 04 = 98.41. 
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FIG. 6: (NPH) The 68% and 95% confidence levels of Q,^h 2 versus 0,2- The left shows the joint likelihood after marginalizing 
over the frequency 04 for the high frequency grid up to 04 < 700 (63 frequencies). If one includes all frequencies in the 
marginalization the constraint on a 2 disappears (right). 



ics of the code. 

In the low frequency regime we obtain two frequencies 
that give an improvement of A\ 2 ~ 6 compared to no 
oscillations. In the next section we will discuss the results 
from an MCMC run for one of these frequencies, with 
a varying phase and amplitude. For the NPH model, 
the results can be summarized as follows: for frequencies 
1 < a 4 < 200, a 2 < 0.13 at 68% and a 2 < 0.21 at 95% 
confidence. For frequencies up to a 4 = 700, a 2 < 0.29 at 
68% and a 2 < 0.39 at 95% while for higher frequencies 



the amplitude is not constrained within the bounds of 
the grid. The most likely value within the low frequency 
grid is ci4 = 46.5 with an amplitude of a 2 ~ 0.056. 

For the BEFT model, there is a constraint set by BEFT 
validity on the value of Aq. A BEFT approach to initial 
state modifications is only valid if the physical momen- 
tum is smaller than the cut-off scale at the boundary, 
i.e., k/ao < M. The maximum value of k is set by the 
smallest observable scale in the CMB, k max = O(0. 1) and 
therefore we deduce Aq < 10 from eq. ([!]). The number 
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FIG. 7: The primordial spectrum from the NPH model with 
a 4 ~ 46, a 2 = 0.14, a 3 = -0.04 and <f> = 1.5. 




FIG. 8: Best fit for the NPH model from the MCMC chain, 
with a fixed frequency <X4 ~ 46 . The improvement of the 
fit compared to no oscillations is A^ 2 ~ 12. From the grid 
we found that an equally good fit is for a frequency 04 ~ 
98, about double this frequency. We find that the best fit 
amplitude is rather large, 0,2 ~ 0.14, with the 68% level still 



allowing zero amplitude (fig. 11 1. The amplitude is similar to 



the best fit amplitude of the axion monodromy model derived 
in [12]. However, this best fit was at a translated frequency 
of about 150. 



of samples in A is set to 120 equidistant values between 
1 and 10. For BEFT the frequency is not constrained 
by slow-roll and is proportional to M/H; the scale of 
new physics divided by the Hubble scale. Consequently 
the effective frequency can be quite high. We sample 700 
logarithmically spaced steps between 10 3 and 10 4 making 
up a total of 1,344,000 grid points. 

The likelihood confidence contours for the BEFT 



model are different from those of the NPH model (fig. 10 1 



Most importantly, the resulting contour does not have a 
vanishing amplitude, which was assumed as a prior based 
on the theoretical form of A . The best fit point in the 
grid is given by A\ = 7708 with A = 10 corresponding 



FIG. 9: The best-fit angular power spectrum derived from 
the grid for the BEFT model with A = 10 and Ai = 7708, 
resulting in A^ 2 — 13 relative to no oscillations. 



to an improvement of A% 2 ~ 13. 



V. MCMC AND MODEL CONSTRAINTS 

We analyzed the WMAP data using Monte Carlo 
Markov sampling with a fixed frequency derived from 
the grids in the previous section for both BEFT and 
NPH model. We set a Gelman and Rubin criterion of 
R— 1 < 0.01 [37]. The Gelman- Rubin diagnostic R relies 
on parallel chains to test whether they all converge to the 
same posterior distribution by considering the variance of 
the parameters in each chain compared to the variance 
of the same parameters over all parallel chains. Conver- 
gence is diagnosed once the chains have 'forgotten' their 
initial values, and the output from all chains has become 
indistinguishable (R — 1 = 0). For the NPH model, we 
ran 8 parallel MCMC with 04 ~ 45.9, close to the best 
fit point in the low frequency grid. In addition to 02 we 
allowed both 0,3 and £ to vary. We assumed flat priors 
with -0.6 < a 2 < 0.6, -0.5 < a 3 < 0.5 and < C < n. 
We did not reach our desired GR criterion with the first 
chains, which is mostly entirely due to the slow conver- 
gence of the amplitude P*. After we ran 8 chains, we 
derived a proposal covariance matrix. We used this co- 
variance matrix to speed-up the convergence. Obviously, 
if the proposal is derived from chains that have not con- 
verged sufficiently, we will not recover the true posterior 
distribution. As such, using a proposal covariance matrix 
is not without risks. With the addition of a covariance 
proposal we easily obtained R — 1 < 0.003 over 4 chains. 
The result for some of the marginalized and joint likeli- 
hoods is shown in fig. [TT] for the NPH model. Note the 
high correlation between a2, the amplitude of the oscil- 
lation, and CZ3, the scale dependence of the oscillatory 
correction. We found no evidence for a strong correla- 
tion between a 2 and 03 and the other parameters, which 
proves that a grid sampling is a very good first estimate 




10000 



FIG. 10: (BEFT) The 68% and 95% confidence levels for amplitude Ao versus the frequency Ai marginalized over the baryon 
density. Again, there are many local maxima in the likelihood. The cut-off in amplitude Ao is set by theoretical constraints. If 
we would allow larger amplitudes (Ao > 10), the likelihood contours would likely move towards higher A. The most likely grid 
point is (Ai = 7708, Ai = 10). 



of the best fit values of both parameters. We found that 
£ is weakly correlated with both the amplitude P* as 
well as the dark matter energy content of the Universe, 
explaining the relatively slow convergence of these distri- 
butions. 

From the MCMC for the NPH scenario we deduce 
marginalized best fit values a 2 — 0.15±0.17, 04 = 
0.04 ± 0.04 and f = 1.6 ± 0.7 with a fixed frequency of 
124 = 45.9 and a best fit point with A\ 2 — 12. Therefore 
we conclude that at the 1 sigma level, 02 is consistent 
with zero. 

To relate these to constraints on the NPH model one 
must investigate the relation between 02, 0,3 and 04 as 
derived in |16| . The exact form of these parameters is 
presented in eq. (31) of their paper 4 We derive a relation 
between 02, 03 and 04 



minus sign. We find M/H > 10 3 and 



«2 
«3 



1 M 

(Z4 H 



[5/2 + In M/H] 



-1 



(4) 



Here we took the lower limit A = M, where we consider 
M/H > 10 2 . The observational limits were derived for 
a 4 ~ 46, which puts a bound on M/H > 10 3 for e < 
0.01 (slow-roll). This results in a theoretical constraint 
0-2/0-3 < —2.3. We can derive a similar constraint for the 
upper limit A = \(H + M 2 /H), or A/H ~ l/2(M/H) 2 . 
It turns out that 04 in this limit can only have a negative 
sign. For a sine, this means the amplitude picks up a 



4 We derived these constraints in the assumption ei/e2 
H„ ~ H due to weak scale dependence 1161 . 



1 and 



02 ~ i^(2-lnM/2tf) 

03 04 H 



5/2 + In -(M/H) 2 



(5) 



which together results in 02/03 > 5.86. 

In fig. [12] we show the joint likelihood together with 
the theoretically allowed values of <Z2 and 03 for a given 
frequency 04. If this frequency is a valid signal, a minimal 
improvement of the confidence contours could exclude 
this oscillation to be due to a NPH altogether. 

For the BEFT model, the bound on the frequency is 
not that stringent and we simply chose the best fit fre- 
quency of Ai = 7708 obtained from the grid. Note that 
this high frequency requires us to compute all I in order 
to resolve the primordial oscillation. This significantly in- 
creases computation time for each power spectrum. From 
the MCMC we would like to derive proper distributions, 
which could be an issue if we put priors on Ao that ex- 
clude values that have non-zero probability. Therefore 
we allow < A < 25 and — tt < <f> < tt. 

The 8 parallel MCMC for the BEFT model resulted 
in an improvement of A\ 2 ~ 16, with all parameters 
satisfying the Gelman-Rubin diagnostic, except Ao and 
</>. We ran 6 additional parallel chains with an estimated 
covariance matrix from the first 8 chains. The Gelman- 
Rubin diagnostic was reduced for both Ao and <J) to R — 
1 < 0.01. 

Again we find no indication that Ao and <fi are strongly 
correlated with any of the other parameters. In fact, Ao 
almost seems independent of the other parameters. This 
is reflected in the fast convergence of the distribution of 
this parameter, which reached a R—l < 0.2 after running 
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FIG. 11: Marginalized and joint likelihoods for several parameters in the new NPH scenario. The parameters ai and az show 
strong correlation. 



the first 8 chains (which is remarkable given the length 
of the chains). Again, the phase 4> is weakly correlated 
with all the energy densities, leading to a relatively slow 
convergence of these distributions, even after estimating 
the covariance matrix. This correlation however is even 
less profound here than in the NPH model. These pa- 
rameters reach R — 1 < 0.02 within the limited running 



time (6 chains, each 120 hours on 8 core 2.2 Ghz CPU's, 
resulting in chains of approximately 45000 samples) of 
the chains. 

For a comparison with theoretical bounds, we derive a 
relation between Aq and A\ as A\ ~ Aq x (M /H) (see 
[H] table 1). For a fixed frequency of A\ ~ 7708 and 
in the assumption the true distribution is a Gaussian 
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FIG. 12: Joint likelihood for ai and 03 together with the the- 
oretical constraints relating a-z to 03 showing 68% and 95% 
confidence levels. The region (shaded, red) allowed within 
the NPH model constitutes a small part of the observation- 
ally allowed region. It shows that a marginal improvement 
of the confidence contours could exclude NPH as the source 
of this oscillation. In addition, both parameters are perfectly 
consistent with zero. 

we derive 7.0 < A < 15.2 at 68% confidence, which 
is almost 3 sigma away from zero. From this value we 
derive that 500 < M/H < 10 3 at 68% confidence which 
lies within theoretical predictions. Based on the best fit 
amplitude A ~ 13 and the average likelihood we obtain 
a slightly different constraint of consequently we obtain 
450 < M/H < 850 at 68% confidence. We argued that 
in principle Aq is constrained to be smaller than 10 in 
the BEFT model. We set the coupling constant (3 — 1 
but in fact j3 is 0(1) which allows for small increase of 
the total amplitude, /3Aq, possibly matching the best fit 
amplitude. Unlike the best-fit oscillation from the NPH 
model, this correction is inconsistent with zero at the 3 
sigma level. 



VI. EFFECTS ON OTHER PARAMETERS 

Although there is little to no correlation between most 
of the parameters from the oscillating component and 
the other parameters within ACDM, the cross correla- 
tion between ACDM parameters is affected differently 
for BEFT and NPH (the best fit values however are simi- 
lar). For example, although confidence the levels of n s are 
hardly affected by the presence of oscillations, in the NPH 
model, weak correlation between the phase and the en- 
ergy densities causes the probability distribution in Q dm 
to broaden, i.e., the correlation between the phase and 
the energy densities in the NPH scenario causes the un- 
certainty in those parameters to increase. This effect is 
shown in fig. |13| 



VII. CONCLUSIONS 

We used the latest WMAP data to constrain oscilla- 
tions on top of an almost scale invariant primordial power 
spectrum. We argued that the primary difficulty in con- 
straining these models with the data is related to the 
irregular likelihood function of the frequency of these os- 
cillations. There are many equally good fits and these 
fits are at discrete frequencies. For all other parameters 
one tries to constrain in a ACDM model of the Universe, 
such degeneracy usually does not exist and it is quite 
sufficient to scan parameter space using an MCMC. In 
order to avoid the random jumps from one maximum in 
the likelihood to another, it is preferable to keep the fre- 
quency fixed. However, the frequency is on of the key pa- 
rameters being constrained, and it would make no sense 
to fix it prior to analyzing the data. In order to estimate 
the frequency before the MCMC we applied a grid search 
for the best fit. This grid quickly becomes incalculable 
for a large number of parameters (which is the primary 
reason to run an MCMC instead) and we only varied the 
parameters characterizing the modification of the power 
spectrum as well as the baryon density, which previously 
had been identified to be correlated with the parameters 
of the oscillatory correction. Once the best fit had been 
established we performed an MCMC with the frequency 
fixed to its best fit value. 

We found that the addition of oscillations on top of the 
smooth power spectrum can improve the overall fit. An 
improvement up to Ax 2 ~ 16 was found once we applied 
an MCMC with a predetermined best fit value of the fre- 
quency, which resulted in a best fit amplitude almost 3 
sigma away from zero. We did not find significant cor- 
relation between the oscillatory parameters besides the 
frequency and the parameters in ACDM, confirming that 
grid sampling is a fairly robust first estimate of the best- 
fit frequency. 

In the NPH scenario very high frequencies should lead 
to unobservable amplitudes. For completeness we an- 
alyzed such high frequencies. We found that for some 
specific frequencies we could obtain an improvement in 
the goodness of fit up to A% 2 ~ 12 compared to no os- 
cillations. Analyzing the angular power spectrum with 
these extreme primordial frequencies revealed that these 
frequencies are in fact unresolved at large angular scales 
(small I), where the wavelength of the oscillating primor- 
dial power spectrum is smaller than the angular distance 
between subsequent Vs. This leads to small glitches in 
the slope of the first peak of the angular power spec- 
trum, which fit some of the observed outliers in the slope 
of the first peak in the observed power spectrum. As 
such, high frequencies, although unresolved, could ac- 
count for some of the large scale effects we observe in the 
data, and improve the overall fit. We conclude, however, 
that these oscillations cannot be caused by NPH modi- 
fications and must be due to a different model, possibly 
an axion-monodromy inflation type model |12j . 
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FIG. 13: The 68% confidence contours between the phase £ (NPH, green) and cj> (BEFT, red) and the dark matter energy density 
and the amplitude. For sake of comparison we shifted the central value of <f> upwards, while contour levels were maintained. 
Both examples show how the correlation between £ and P„ and Q.dm affects the uncertainty in these parameters. It does seem 
that overall uncertainty is conserved, as the error in <f> is larger than the error in f. 



Although in the grid for the low frequency domain 1 < 
£14 < 200 of the NPH model we could only found one 
frequency with an improvement of A\ 2 — 6, we ran an 
MCMC near this frequency to obtain an improvement 
of Ax 2 — 12 corresponding to a best-fit amplitude of 
a-i = 0.16. This is close to the best- fit amplitude found 
in [12 with an improvement of A\ 2 ^11 (WMAP-5). 
Their improvement appeared at a (translated) frequency 
of 04 ~ 150. It is interesting to note that there appear to 
be subsequent improvements at low frequencies close to 
50, 100 (WMAP-7) and 150 (WMAP-5) which could be 
an additional hint we might be looking at an oscillation 
as opposed to noise, because the improvement appears 
at equidistant intervals in frequency. 

Using an MCMC we have been able to put constraints 
on several primordial parameters. In particular, we de- 
rived constraints on the amplitude a 2 , the tilt 03 and the 
phase C of the oscillatory correction in the NPH model. 
We used these observational constraints to test the NPH 
model by implementing the relation between various pa- 
rameters as predicted by the NPH model. This shows 
that the signal found at this frequency could originate 
from NPH modifications to the primordial power spec- 
trum. However, an improvement of the confidence levels 
would exclude this possibility. 

For the BEFT model we ran one grid and we found 
that even before varying all parameters in an MCMC we 
could achieve an improvement of Ax 2 — 13. This im- 
provement corresponded to a frequency A\ = 7708 and 
an amplitude Ao = 10. We ran an MCMC around this 
frequency and found that the theoretical bound on A 
does not probe the observationally best fit point. In or- 
der to recover the full distribution of Aq we set an upper 



limit on A < 25 when running the MCMC. We found 
a best fit improvement of A\ 2 ~ 16 compared to no os- 
cillations with a phase cp — —0.21 and an amplitude of 
02 = 12.25. This best fit value of Ao is theoretically 
not expected by a BEFT modification of the primordial 
power spectrum. The amplitude of the oscillatory correc- 
tion is also proportional to a coupling constant /3 which 
we set to 1. Once we relax this assumption, larger am- 
plitudes are theoretically possible through f3A a . These 
results enabled us to put a constraint on the ratio M/H 
of 5 x 10 2 < M/H < 10 3 , which is a very realistic possi- 
bility of this parameter ratio. 

It is hard to assess the significance of the improved fits. 
We could think of two possibilities: investigate simulated 
CMB data without oscillations with realistic beam and 
noise, and investigate simulated CMB data with oscilla- 
tions and realistic beam and noise. For example, if we can 
find oscillations in simulated data without oscillations, 
this would suggest oscillations can well be mimicked by 
noise. Complement arily, generating realistic CMB data 
including primordial oscillations, can we recover these os- 
cillations by analyzing the data? In particular, it would 
be interesting to determine whether there exists a thresh- 
old amplitude for oscillations to be recovered. Preferably 
we should generate a large number of maps in order to 
see if we can assess a probability of recovering oscillations 
from the data. We will report our findings in future work. 

For now, we conclude that a primordial power spec- 
trum with no oscillations is consistent with the data for 
most amplitude and frequency of the primordial signal. 
For some frequencies, a non-zero amplitude of the os- 
cillatory corrections to the power spectrum seems to be 
preferred by the data. These signals can be investigated, 
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and used to constrain primordial parameter space, possi- 
bly signifying some of the detailed physics driving infla- 
tion. 
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